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Abstract 

We review recent results from computer simulation studies of polymer glasses, 
from chain dynamics around glass transition temperature Tg to the mechanical be- 
haviour below Tg. These results clearly show that modem computer simulations 
are able to address and give clear answers to some important issues in the field, in 
spite of the obvious limitations in terms of length and time scales. In the present 
review we discuss the cooling rate effects, and dynamic slowing down of different 
relaxation processes when approaching Tg for both model and chemistry-specific 
polymer glasses. The impact of geometric confinement on the glass transition is 
discussed in detail. We also show that computer simulations are very useful tools 
to study structure and mechanical response of glassy polymers. The influence of 
large deformations on mechanical behaviour of polymer glasses in general, and 
strain hardening effect in particular are reviewed. Finally, we suggest some direc- 
tions for future research, which we believe will be soon within the capabilities of 
state of the art computer simulations, and correspond to problems of fundamental 
interest. 

1 Introduction 

A polymer is a macromolecular chain resulting from the connection of a large number 
of monomeric units. The number of monomers (the chain length A^) typically ranges 
between 10^ and 10^ in experimental studies of polymer melts [JJ. Recently, this range 
of chain lengths has also become accessible in simulations IIH O . Such a long chain 
has an open structure which is strongly pervaded by other chains in the melt [4 |. The 
strong interpenetration of the chains has important consequences for the properties 
of the melt. For instance, intrachain excluded volume interactions, which swell the 
polymer in dilute solution ||T]|4l, are almost screened so that a chain behaves on large 
length scales approximately as a random coil [3J. Furthermore, chain interpenetration 
also impacts the polymer dynamics by creating a temporary network of entanglements. 
Entanglements strongly slow down the chain relaxation and make the melt viscoelastic 
already at high temperature H] |6l . 
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Figure 1 : Volume-temperature diagram of a coarse-grained model for poly(vinyl al- 
cohol) [7,i8J. The volume per monomer (v) and the temperature (T) are given in 
Lennard- Jones units. In the liquid phase the chains have random-coil-like configura- 
tions and the structure of the melt is amorphous. The amorphous structure is preserved 
on rapid cooling and the melt eventually undergoes a glass transition at low T. For 
slow cooling the melt transforms into a semicrystalline material in which sections of 
folded chains order in lamellar sheets that coexist with amorphous regions. When the 
crystal is slowly heated up (dashed line), it melts at a higher T than the temperature 
where crystallization occurs. This hysteresis is typical of first-order phase transitions. 

On cooling toward low temperature the polymer melt eventually transforms into a 
solid. Polymeric solids are either semicrystalline or glassy (Fig.[T]) |9|. In the semicrys- 
talline state amorphous regions are intercalated between crystalline lamellar sheets. 
The sheets consist of chains which are folded back on themselves so that chain sec- 
tions align parallel to each other |i9jl10J. The ability to form crystals crucially depends 
on the polymer microstructure. Only chains with regular configurations, e.g. isotatic or 
syndiotatic orientations of the sidegroups or chains without sidegroups, can fold into 
crystalline lamellae. However, even in these favorable cases full crystallization is hard 
to achieve ifTOl . 

This intrinsic difficulty to crystallize favors glass formation ifTTlfTSl fTTl. Polymer 
melts either can be easily supercooled (Fig.lTJ or, due to the irregular chain structure, 
do not crystallize at all. An example for the former case is bisphenol-A polycarbon- 
ate (PC) 1 14|, while examples for the latter involve homopolymers with atactic (bulky) 
sidegroups, such as atactic polystyrene (PS), or random copolymers, such as cis-trans- 
1 ,4-polybutadiene (PBD). These polymeric glass formers share with other (intermedi- 
ate and fragile) glass-forming liquids the key characteristic feature of the glass transi- 
tion; that is, little change of the amorphous structure, but a huge non-Arrhenius-like 
slowing down of the dynamics on cooling toward the glass transition temperature Tg 
lfT2l [131 [T5l . Understanding the molecular origin of this disproportionate behaviour 
represents a great scientific challenge 1121 [131 [TSl [T6l flTi . In addition to this funda- 
mental interest in the study of the glass transition, soUd polymers are also integral com- 
ponents in many modern applications [18|. For instance, glassy polymeric materials 
are appreciated because of their spectacular mechanical properties |[T9ll20ll21J . Instead 
of failing abruptly when subject to strong deformations, some polymers, such as poly- 
carbonate, may harden for large strains, leading to a tough mechanical response. A 
microscopic understanding, elucidating the structure-property relationship, of this be- 
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haviour is still elusive. A further example is provided by thin polymer films, extensively 
used in technological applications as protective coatings, optical coatings, adhesives, 
etc. In the ongoing quest for progressively smaller structures and devices these films 
may attain nanoscopic dimensions, where deviations from the bulk behaviour should 
be expected (2T\. Indeed, many recent studies suggest that the of thin films is shifted 
by spatial confinement, a striking observation which is not well understood yet 1 23 24| . 

This wide array of challenging questions — from the glass transition to the impact 
of external stimuli — ^bestows the theoretical understanding of glassy polymers with 
particular significance. Molecular simulations can contribute to this research. Over 
the past two decades, computational models and methods have been developed for 
simulating these glassy systems. The progress made in the field has regularly been 
the subject of topical reviews. For instance, detailed reports of chemically realistic 
modeling approaches may be found in Refs. Il25ll26ll27l . while work on coarse-grained 
models is reviewed in Refs. 1281 l29l [30l [311 [32l |27]| . The purpose of the present article 
is to give a brief account of recent work. 

We begin our survey by a short introduction to the modeling of polymers in sim- 
ulations (section |2|. The glass transition temperature naturally splits the following 
discussion into two parts, a part devoted to the properties of the model glass formers 
above Tg (section[3]) and a part addressing sub-Tg phenomena (section|4|i. In section|3] 
our discussion is mainly concerned with the dynamics of weakly supercooled polymer 
melts and the impact of spatial confinement on their behaviour. Section |4] focuses on 
the response of polymeric glasses to both weak and large external deformation, where 
respectively the linear mechanical behaviour of the glassy melt and its approach to 
material failure will be explored. The article concludes with a short summary of the 
presented results and an outlook on possible future research directions. 

2 Computer simulations: models and computational as- 
pects 

Molecular simulations of glass-forming polymers face the important challenge that 
the structure and dynamics of these systems are governed by a large spread of length 
and time scales |33, 34|. The relevant length scales extend from the atomic diameter 

10^^" m) to the chain dimension 10^^ m for ^ 10^). The spread of time 
scales is even larger; it ranges from bond vibrations 10^^^ s) to the slow structural 
relaxation close to Tg 10^ s). 

From the inspection of these scales, it is clear that simulation approaches have some 
limitations and also require simplifications. An obvious limitation are the accessible 
time scales. With the currently available computer power the longest time that can 
be reached in molecular dynamics (MD) simulations is roughly a few /is. This time 
scale is 8 orders of magnitude smaller than the relaxation time at the experimental Tg, 
implying that the simulated Tg is shifted by about 25 degrees to higher temperature (T) 
relative to experiments (where we used the rule of thumb: 1 decade in time S 3 K). The 
simplifications that are currently necessary concern the simulation models which are 
obtained by some kind of coarse-graining procedure, designed to eliminate (some) fast 
degrees of freedom by incorporating them in effective potentials. This model-building 
step can take different levels of complexity 1*3511361121137 1, but roughly speaking, there 
are two families, atomistic models and coarse-grained models. In the following, we 
briefly introduce these models and discuss their strengths and weaknesses as we go. 
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Before doing so, however, we want to point out some general appealing features of 
the simulations to balance the reservations expressed above. The simulations offer full 
control over the perfectly defined system under study. For instance, it is possible to vary 
in a systematic manner parameters, such as chain stiffness, chain length, etc., while 
keeping all other defining properties of the model, so as to single out the impact of each 
parameter on the properties of the system. Furthermore, the simulations enable one to 
study local properties or to explore correlations which are hard (or even impossible) to 
access in experiments. This provides microscopic insight and a means to test theoretical 
concepts. Examples from the recent literature involve investigations of the potential 
energy landscape of glass-forming materials Il38ll39ll40l . of single-molecule diffusion 
in polymeric matrices 14111421 or of fracture in glassy polymers IZTl 

Atomistic models. Atomistic models replace electronic degrees of freedom by force 
fields, i.e., empirical potentials for bond-length, bond-angle, torsional-angle, and non- 
bonded interactions, whose parameters are determined from quantum-chemical cal- 
culations and experiments ll26l [35l . During the past decades, force fields have been 
developed for both explicit atom (EA) models, treating every atom present in nature 
as a separate interaction site, and united atom (UA) models which represent a small 
number of real atoms (e.g., CH2 or CH3) by a single site. The reduction of the num- 
ber of interaction sites in the UA model has the computational advantage of allowing 
longer simulation times. With a time step of ^ 10^'^ s a few thousand united atoms 
can be simulated over several 100 ns, about an order of magnitude longer than an EA 
simulation of comparable system size. 

Both EA and UA models were employed to study the properties of glassy poly- 
mers (see e.g. Il25l l26l l27l for reviews). Recent examples include polyisoprene (EA 
model (4T. '441), atactic PS (UA model (AF, '46', '4T, '48', '491), bisphenol-A PC (UA 
model 148, 49J) 1,4-PBD (UA and EA models |26, 50, 5L, 52, 51 |55] |56j), and 
poly (ethylene oxide) or atactic poly (propylene oxide) (EA models If57l ). Certainly, the 
main strength of these modeling efforts is that the simulation results allow for a direct 
comparison with experiments. In some cases, such a comparison is possible with com- 
mercial software packages ll56l . while a careful fine-tuning of the force field is often 
required fSO^/S?!. Experience from those modeling approaches — not only for poly- 
mers but also e.g. for amorphous Si02 L3_l , 5S| — suggests that the design of a chem- 
ically realistic model, aiming at a parameter-free comparison between simulation and 
experiment, should involve information about both structural and dynamic properties 
Il26i. 

Generic models. Atomistic simulations are ideally suited to study specific polymers, 
including their glass transition and the properties of the glassy state. On the other 
hand, the strong increase of the relaxation time, which eventually leads to vitrifica- 
tion on cooling through Tg, is common to all glass-forming polymers, irrespective of 
their chemical composition and architecture. This universal aspect suggests to em- 
ploy simplified simulation models which only retain generic features of a polymer 
Various such generic models — on a lattice [ 3T1 [30l l28l l28l or in the spatial contin- 
uum 1 3TI [30I [32I — ^have been studied for glass-forming polymers. In the continuum, 
the simplest model are highly flexible bead-spring models, where spherical monomers 
are tethered together by springs and have nonbonded, Lennard-Jones (LJ) interactions 
IISTl l30l [32I . More chemical realism can be introduced by making the chains semi- 
flexible, through the addition of a bond-angle potential ||59l and, possibly also, of a 
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torsional potential ll^ ISTl l62l . Due to their simplicity, these generic models allow for 
an efficient simulation. This is important if one wants to vary the cooling rate |63| or 
the strain rate |27 1 over decades, to explore systematically the impact of model parame- 
ters, such as chain length or chain flexibility Il59ll60l 16111621 . or to obtain good statistics 
for comparison with theory [64\. Computational expedience is also important for ex- 
ploratory studies of more complex systems, such as inhomogeneous polymer systems 
(e.g., polymer films, nanocomposites, semi-crystalline polymers) or glassy polymer 
mixtures (e.g., dynamically asymmetric polymer-polymer mixtures or polymer-solvent 
systems). Some examples will be discussed in Sects. [3]and|4] 

Hierarchical models. For generic models, the gain in the accessible length and time 
scales is obtained at the expense of a loss of correlation to the atomistic conforma- 
tion of the polymer For many problems in materials research, this loss is undesirable. 
Therefore, much research efforts currently goes into the development of hierarchical 
approaches consisting of interconnected levels of modeling (atomistic, generic, macro- 
scopic) lIJl [37] |65] 166). The idea is that each level treats phenomena on its specific 
length and time scales and then passes on the results as input to the next, more coarse- 
grained level, until the desired materials properties can be predicted. Such multiscale 
simulation methods represent a powerful approach whose potential for the modeling of 
glassy polymers is beginning to be explored (see e.g. 167)1 ). 

Remarks on simulation methods. Sections[3]and|4]will present results from molec- 
ular dynamics (MD) simulations, a numerical method to integrate the classical equa- 
tions of motion for a many -body system in a given thermodynamic ensemble I68l l69]| . 
Therefore, MD is the natural simulation technique to address dynamical problems, 
such as the glass transition. However, the realistic MD dynamics carries an obvious 
price: the equilibration time for the system under consideration — for instance, for a 
long-chain polymer melt close to its Tg — will exceed the maximum time of a few mi- 
croseconds one is currently able to simulate. 

Here Monte Carlo (MC) techniques may provide a promising avenue because of 
large freedom to design MC moves 1681 iTOl . The hope is to find an efficient algo- 
rithm allowing one to decorrelate the configurations of glassy polymer melts rapidly 
For long-chain polymer melts, this demand on the algorithm implies, already at high 
T, that the MC move should be nonlocal, i.e., it should modify the chain conformation 
at large scales, and it should not require empty space because the melt is a dense liquid. 
A promising algorithm satisfying these requirements employs double-bridging moves 
which alter the connectivity between two neighbouring chains while preserving the 
monodispersity of the chains (for recent reviews see e.g. 13] |2])- Such a connectivity- 
altering move drastically changes the conformation of the two chains involved and thus 
relaxes the length scales on the order of the chain dimension efficiently. However, it 
does not alter the local packing of the monomers. An inherent hazard of the algorithm 
therefore is that, if the move is attempted repeatedly, a successful double-bridging event 
is likely to annihilate one of its predecessors by performing the transition between two 
chains in the reverse direction. To avoid this inefficiency the nonlocal chain updating 
should be complemented by a move which efficiently mixes up the local structure of the 
melt. At low T, efficient relaxation of the liquid structure calls for a method which al- 
leviates the glassy slowing down in general. Thus, any algorithm achieving this aim in 
nonpolymeric liquids should also accelerate the equilibration of glassy polymer melts, 
provided that it can be generalized to respect chain connectivity. At present, no tech- 
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nique has been established to solve this problem. However, possible candidates could 
be parallel tempering 117111721 1731. Wang-Landau sampling |70| or variants thereof |74|, 
or transition path sampUng methods 1,75 J . 



3 Glass transition and properties of the supercooled poly- 
mer liquid 

3.1 Bulk properties 

3.1.1 Glass transition temperature 

In polymer melts, the transition from the glass to the liquid is accompanied by a strong 
decrease in the shear modulus, typically of three to four orders of magnitude [21]. It 
is thus clear that the glass transition temperature is a characteristic of high engineering 
relevance. Classical methods for the experimental determination of Tg are calorime- 
try and dilatometry. Both methods hint at the kinetic features of the glass transition. 
The transition occurs when the relaxation time for volume recovery (dilatometry) or 
enthalpy recovery (calorimetry) becomes longer than the time scale of the experiment 
(i.e., than the cooling or heating rate) IfTTl [121 [T3l . Upon cooling the polymer liquid 
falls out of equilibrium close to Tg and freezes below Tg in a glassy state, the properties 
of which depend on the details of the cooling process and tend to age physically during 
further isothermal equilibration |[lTl[T2l[T3l . The glass transition temperature is located 
in the T interval where the polymer melt smoothly evolves from the liquid to the solid 
state, and can be defined operationally through some prescription, for instance, as the 
intersection of straight-line extrapolations from the glassy and liquid branches of the 
volume-temperature curve [flTllTSl . 

Similar extrapolation procedures are also applied in simulation studies Il26ll63l[30l 
|42l|48||77||28||76||79]|80l. The resulting Tg values have dependences comparable to ex- 
perimental ones despite the much larger cooling rate employed — typically lO'^K/min 
in simulations and lOK/min in experiments. For instance, Soldera and Metatla find a 
linear relationship between numerical and experimental Tg values from atomistic simu- 
lations of various vinylic polymers |77|. Figure |2]reveals that Tg decreases nonlinearly 
with the logarithm of the cooling rate (Fj-) Ii26..63.,30ii47j . 

r„(r7') = 7;'' — - — -, (i) 

' s in(Arr) 

in accordance with experimental observation lISTl . Also in agreement with experiment 
[[nilHllEllp 851, Tg increases with chain rigidity and chain length |60l|47l[86l. The 
chain length dependence can be fitted to the empirical Fox-Flory equation 

TsiN)=Tg-^, (2) 

which usually describes experimental data well (if the molecular weight is not too small 
lim ). although other forms have recently been discussed in the literature |[84l[85ll and 
can be rationalized theoretically liSTl . 

3.1.2 Dynamics of the supercooled melt 

The precursor of the glass transition is the strong slowing down of structural relaxation 
processes on approach to the transition from the liquid. This dynamical feature is 
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Figure 2: Panel (a): Glass transition temperature (Tg) versus cooling rate (Tj) for 
an atomistic model of atactic PS |47|. The simulation box contains 8 chains with 
= 80 monomers each. The solid line shows a fit to Eq.[T]with To = 371 K, B = 110 
K, and A — 0.23 ps/K. Panel (b): Tg versus chain length (A^) for a fully flexible and 
a semiflexible (with angular potential) bead-spring model [ f76l . The simulation box 
contains at least 192 chains and the cooling rate is Fr = 2 x 10^^. The solid lines are 
results of a fit to Eq.|2]with Tg = 0.432, K = 0.145 (flexible model) and T^ = 0.525, 
K — 0.264 (semiflexible model). All data are given in Lennard-Jones units for the 
bead-spring models. 
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t[ps] 

Figure 3: Mean-square displacement (MSD) guit) versus time t for a united atom 
model of 1,4-PBD (results adapted from |53|, with permission). The MSD is averaged 
over all united atoms in the melt containing 40 chains of 30 repeats units. The temper- 
atures shown are (from left to right): T = 400 K, 323 K, 273 K, 240 K, 225 K, 213 
K, and 198 K {Tc ~ 214 K). The horizontal dashed lines indicate the radius of gyration 
/?2~218 A2, the average Lennard-Jones diameter of the united atoms <7 w 3.8 A, and 
an estimate for the Lindemann localization length r^c ~ 0.45 A ||43]| . The dotted lines 
for T = 400 K and 213 K show the power law ^ f"-^', characteristic of Rouse-like 
motion. 
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a hallmark of strongly interacting disordered matter, including dense colloids ||88]| . 
granular materials [89|, and supercooled liquids 1 12, 16 1. In simulated polymer melts, 
it is observed in all dynamic correlation functions, e.g., in dynamic structure factors, 
conformational correlation functions, or dielectric relaxation |26|. 

As an example. Fig. |3] shows the mean-square displacement (MSD) guit), aver- 
aged over all monomers of a chain, for a chemically realistic model of 1,4-PBD Il53l . 
At high temperature, the MSD directly crosses over from ballistic motion t^) at 
short times to subdiffusive motion (~' f^" with xq ~ 0.61) at intermediate times where 
the MSD is bound between the monomer size (a) and the end-to-end distance of a 
chain. This subdiffusive motion does not depend on the strength of the torsional barrier 
(cf. Fig.m, is present even if the torsional potential is absent |52|, and thus reflects the 
universal Rouse-like dynamics of nonentangled chains in a polymer melt |4|. The term 
"Rouse-like" stresses the fact that simulations of nonentangled chains |26, 64 901I9T1 
find both accord with Rouse predictions — e.g., the Rouse modes are (nearly) orthog- 
onal for all t — and deviations from them — e.g., the Rouse modes are stretched with 
stretching exponents depending on the mode index and the MSD of the chain's center 
of mass increases sublinearly for intermediate times. The origin of these deviations is 
not fully understood. However, for short chains the relaxation time of a chain is not 
well separated from the (local) a-relaxation so that finite-A^ corrections to the Rouse 
behaviour must be expected |64|. Moreover, even for long chains intermolecular in- 
teractions between the polymers are not completely screened on mesoscopic length 
scales, which may cause subdiffusive center-of-mass motion ||92l or deviations from 
reptation theory |6|, and lead to corrections to chain ideality If3l l93l . 

On cooling towards To the Rouse-like motion shifts to progressively longer times 
due to the appearance of a plateau regime. In this regime, the MSD increases only 
very slowly with time and is of the order of 10% of the monomer diameter, reflecting 
the temporary localization of a monomer in the cold melt. For polymer models with 
intramolecular rotational barriers, such as PBD, this intermittence of large scale motion 
stems from two mechanisms of dynamical arrest Il26l ISTl p5l : monomer caging by 
near neighbours in the dense melt and the slowing down of intrachain conformational 
transitions occurring through (correlated) torsional motion. 

Recent work by Smith and Bedrov |54| suggests that the interplay of both mecha- 
nisms can explain the Johari-Goldstein j3 relaxation in polymers 1 13]. Figure |4]repro- 
duces one of their results. Since experiments reveal that the separation of the a and j3 
processes occurs only on time scales significantly longer than the multiple microsecond 
trajectories generated for the chemically realistic (CR) PBD model. Smith and Bedrov 
employ a PBD model with reduced torsional barriers (LB model) but otherwise iden- 
tical interactions as for the CR model. While structural properties of the PBD melt 
remain unaffected 115211511 . the reduction of the barriers accelerates the dynamics and 
shifts the a-/3 bifurcation into the simulation time window. Figure |4] compares g^it) 
with the decay of the torsional autocorrelation function (TACF). Similar to gu{t), the 
TACF relaxes in two steps. To both steps can be associated relaxation times, and 
Tc which display an Arrhenius (t^) and a non-Arrhenius (la.) increase with decreas- 
ing r, characteristic of j3 and a processes, respectively. At low T, both processes are 
well separated from each other, the beta process corresponding to times for which the 
monomers displace, on average, by about 10% of their diameter This implies that the 
cages imposed by the polymer matrix remain largely intact during the j8 process and 
create a potential energy landscape for the underlying conformational transitions. In- 
deed, detailed analysis reveals that (nearly) all dihedrals visit all torsional states during 
the j3 relaxation. However, these visits do not occur with the equilibrium probabil- 
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Fi gure 4: Monomer MSD gM{t) "ind torsional autocorrelation function (TACF) versus t 
for 1,4-PBD. The data for T = 130 K and 140 K (lines) are obtained from a model with 
reduced torsional barriers (LB), but which is otherwise the same as in Fig. |3] (results 
adapted from [541 , with permission). For comparison the symbols show the MSDs at 
r = 198 K and 213 K for the chemically realistic (CR) model with the full torsional 
potential. These temperatures are approximately at the same distance to Tg (« 170 
K) as for the LB model (Ta w 102 K) |54|. The horizontal dashed lines indicate the 
average Lennard-Jones diameter of the united atoms a w 3.8 A and an estimate for 
the Lindemann localization length rsc ~ 0.45 A ['431. The solid line for T = 140 K 
shows the power law ^ f"^', characteristic of Rouse-like motion. The vertical dotted 
lines represent the /3 -relaxation times (Tp = 716.5 ps for T = 130 K, Tp ~ 319.4 ps for 
T = 140 K). 
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Figure 5: Temperature dependence of the a relaxation time (Ta) for united atom models 
of atactic PS and bisphenol-A PC. Za is normalized by its value at T = 600 K for both 
polymers (PS: t(600K) = 1.4 ps, PC: t(600K) = 1.3 ps) and is obtained from the 
monomer MSD in simulations of 1 chain with = 80 for PS and 64 chains with = 10 
for PC. For both polymers, the increase of la at high T is compatible with the MCT 
prediction Ta - (r - T^)-y (PS: Ti ~ 375 K, 3; PC: Ti ~ 450 K, 7~ 2.2). However, 
the divergence at is avoided and a crossover to an Arrhenius behaviour (PS: ~ 30 
kJ/mol; PC: E^— 16.8 kJ/mol) occurs close to Tg (derived from dilatometry. PS: Tg ~ 
375 K; PC: Tg ~ 433 K). This crossover indicates the transition from the cooperative 
translational dynamics above Tg to activated (/3 process like) hoping below Tg. Results 
adapted from f48l. 

ity because the propensity of a dihedral to return to a preferred conformational state 
increases on cooling, due to the stiffening of the matrix. Equilibrium occupancy of 
conformational states is only achieved for times comparable to the a time scale when 
the "cage effect" of the matrix fully decays (see also |56, 57| for similar results). 

In glass physics, the term "cage effect" is intimately connected to the mode-coupling 
theory (MCT) of the glass transition, which has certainly been one of the most influen- 
tial theoretical approaches in the field during the last twenty years. A comprehensive 
review of its foundations and applications was published recently IfTTl . By correlat- 
ing the equilibrium structure to the dynamics of a glass-forming liquid MCT provides 
a framework for interpreting spatio-temporal correlations measured in experiment or 
simulation. Analytical predictions are derived in the vicinity of an ideal glass transi- 
tion which occurs at a critical temperature and is driven by the mutual blocking of 
a particle and its neighbours ("cage effect"). Extensive tests by experiments and simu- 
lations indicate that, although (the extrapolated) lies above Tg and thus no structural 
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Figure 6: Simulation results for a flexible bead-spring model; R^{r,t) versus r for 
different times from the plateau regime at T = 0.48 (7^ ~ 0.45). Panel (a) shows Ry^ (r, f ) 
for the self-part of the van Hove function and panel (b) for its distinct part. In panel 
(b), the pair-distribution function g{r) is also shown (dotted line; rescaled to fit into 
the figure). The first peak of g{r) reflects the bond length 0.97) and the minimum 
position of the Lennard-Jones potential (« 1.12). In panel (a), r\ (= 0.2323) denotes 
the zero of 7? s ('", f ) and the circles represent a Gaussian approximation which has a zero 
at \/6rsc — 0.2327 and a minimum around -y/TOrsc — 0.3. Here, rsc (— 0.095) is the 
Lindemann localization length. In both panels, the dash-dotted lines correspond to the 
time closest to t' which is the least precise because Rxif,t) is undetermined for t = t' . 
Adapted from ||94l|32l. 
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arrest is observed at Ti (cf. Fig.|5]l, MCT still describes many dynamical features well 
lfT7ll95l . Why this is so, represents a great challenge for the theoretical understanding 
Il88]|96l|97l'98 99 88 100|. 

A key prediction of MCT is that structural relaxation functions should obey a fac- 
torization property in the plateau regime (also called j3 regime in MCT). For the self 
(Gs) and distinct (Gd) parts of the van Hove correlation function I lOU this reads 

Gx(r,f)=iS.W+i/xWG(0 (x = s,d). (3) 

The factorization property refers to the fact that the correction to the nonergodicity 
parameter Fx{r) splits into two factors, of which G{t) depends only on time (and tem- 
perature) and //x(r) only on r Ifr7ll95l . Therefore, the ratio (r' = constant) 11102111031 

Gx(r,r)^Gx(r,0 _ H.jr) 
''''^'''^-Gx(r',0-Gx(r',O"//x(r') 

should be independent of t. For a flexible bead-spring model ll94l Fig.|6]confirms this 
prediction directly from the simulation data (no fit) and additionally reveals the length 
scales pertinent for the dynamics in the plateau regime because distances for which 
i/x('") is zero will not contribute to the relaxation. For the self part of the van Hove 
function the dynamics involves displacements up to the monomer diameter, whereas 
for the distinct part it includes monomers up to about the forth neighbour shell. This 
local character of the relaxation is a direct evidence for the cage effect. 

The agreement between MCT and simulation, demonstrated in Fig.|6] is not limited 
to the flexible bead-spring model. Simulations of semiflexible bead-spring models 
with rotational barriers lf6n[62l and of chemically realistic models, including 1,4 PBD 
Il53ll50ll55]| . PS and PC f4S\, generally find that many features of the spatiotemporal 
relaxation in weakly supercooled polymer melts (i.e., T > Tc) are well described by 
MCT (see however Il57i ). The main difference between these models and the flexible 
bead-spring model is that two coexisting arrest mechanisms — intramolecular barriers 
and monomer caging — determine the structural relaxation of the former, whereas only 
caging operates for the flexible model. This interpretation is a key result of Refs. IISTl 
l52i and of the recent work by Bernabei et al lISTl l62l who also make the interesting 
conjecture that, within MCT, the pertinent theoretical framework for polymer models 
with internal rotational barriers are (so-called) higher-order transition scenarios |17|, 
as it appears to be the case for other systems with distinct arrest mechanisms, such as 
dense colloidal suspensions with short-range attraction I104II1051 or polymer mixtures 
with strong dynamic asymmetry II 1061 11071 [T08 j . 

A distinctive feature of MCT is that the dynamics (beyond the short time regime) 
is fully specified in terms of the liquid structure. This opens the way for an ''ab initio'' 
prediction of the simulated dynamics, solely based on static input obtained from an in- 
dependent simulation of the studied glass former Such an atomistic theory for the slow 
relaxation of (nonentangled) polymer melts has been developed II 10911 and compared to 
simulations of a bead-spring model with = 10 111 101 164| . The comparison gives 
semi-quantitative agreement between simulation and theory (for models with large ro- 
tational barriers further complications might arise |62|). As an example. Fig. |7] shows 
various MSDs, revealing the strengths and weaknesses of the approach |64|. Certainly, 
the agreement for the monomer dynamics is very good. MCT describes the dynamics 
in the plateau regime, the following polymer-specific subdiffusive increase, ^ t'^'^^ 
which is identified as a finite-A^ deviation from Rouse behaviour, and the faster motion 
of the end monomer (gi) relative to the central monomer (gs) of the chain. On the 
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Figure 7: Simulation results for a flexible bead-spring model: Log-log plot of the 
monomer MSD (labeled M, left scale) and the MSD of the chain's center of mass 
(labeled C, right scale) versus Dt with D being the diffusion coefficient of a chain. The 
inset shows the ratio gi{t)/g5{t) (end monomer MSD over middle monomer MSD). 
The circle refer to the MD results at T = 0.47 (Tc ~ 0.45), the solid lines to the MCT 
a master curve, and the dashed lines to the MCT predictions including the MCT j3 
process. The dash-dotted lines indicate the diffusive motion 6Dt and the dotted line 



shows the power law 



,0.63 



of Rouse-like motion. Figure taken from 1641 . 
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other hand, the theory is not so satisfactory for the MSD of the chain's center of mass 
(gc)- Besides underestimating the plateau height it does not reproduce the subdiffusive 
center-of-mass motion ||92l (gc ~ f~'' **) between the plateau and diffusive regimes. 

The representation in Fig. |7] — that is, plotting the data versus Dt with D being the 
chain's diffusion coefficient — facilitates the comparison of the time dependence in the 
P and early a-regimes, but camouflages a systematic deviation between theory and 
simulation. MCT predicts that the a relaxation time Ta{q) has the same T depen- 
dence for all wave vectors q, whereas simulations, not only for flexible bead-spring 
models [32J but also for other glass formers [111 1121 11131 1114! 11151 , find that Ta for 
wave vectors smaller than the position q* of the first maximum of static structure factor 
increases on cooling more weakly than for q^q* ■ This difference also implies a vio- 
lation of the Stokes-Einstein relation 111 16111 ITlll 18L i.e., the product of the diffusion 
coefficient (corresponding to the limit <7 — > 0) and Ta{q*) is not independent of T . 

This decoupling has been interpreted as a signature of increasingly heterogeneous 
dynamics in the liquid nearits glass transition 1 117, 118, 119, 120|). The term "hetero- 
geneous dynamics" means that a glass former near Tg contains subensembles of parti- 
cles with enhanced or reduced mobility relative to the average. To reveal this dynamic 
heterogeneity various methods were deployed 1 1 19 120 j, including filtering techniques 
to track slow or fast particles, analysis of ensemble-averaged three- or four-point cor- 
relation functions 01211 11221 . or scrutiny of the self-part of the van Hove correlation 
function 0111llll2l[TT5l . 

Here we briefly discuss the approach of Refs. 1111111121 which has recently been 
applied to simulation data of a flexible bead-spring model P123V Following 111 1111 12l 
we define the probability distribution P(lnr;f) of the logarithm of monomer displace- 
ments in time t by 

P(lnr;f) =47rr^Gs(r,f) , (5) 

where G^{r,t) is the self-part of the van Hove function, as before. This probability is 
shown in Fig. [8] together with guit) and the non-Gaussian parameter II111II1121 

where r{t) denotes the position of a monomer at time t. In the cold melt there are 
clear deviations from Gaussian behaviour at all but the shortest and longest times. The 
non-Gaussian parameter is positive and has a maximum at time in the late-a regime. 
At this time the distribution P(lnr;f) is very broad, exhibiting small (r > \/6rsc) and 
large (r > 1) displacements. This hints at large, non-Gaussian fluctuations in particle 
mobility when peaks. On cooling toward the broad displacement distribution 
develops into a double-peak structure, indicative of large disparities in particle mobility. 
Apparently, two populations of monomers coexist, "slow" ones which have not moved 
much farther than 10% of their diameter in time f^, and "fast" ones which have left their 
cage and covered a distance of about their diameter or more. This bimodal character 
of the structural relaxation is hard to predict from MCT 111 121 . could be responsible for 
the absence of the divergence of the relaxation time at 7^, and appears to be a general 
feature of materials close to glass or jamming transitions 1,1 15J . 
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Figure 8: Simulation results for a flexible bead-spring model: Probability distribution 
P(lnr;f) of the monomer displacements (left scale) at T = 0.43 (7^ ~ 0.415) for f = 
10, 100, 250, 500, 1000, 1500, 2500, 5000, 10000, and 25000 (from left to right). 
The thick dotted line shows P(lnr;f) at f = fj, the peak time of a2{t). The dotted 
curve for t = 25000 represents the Gaussian approximation for f (lnr;f). The vertical 
dashed line indicates the displacement corresponding to the Lindemann localization 
length v^^sc =0.232. Right scale: Monomer MSDgM(0 and non-Gaussian parameter 
0^(f) at r = 0.43. The vertical dotted line indicates the peak time of oc^(f) and the 
horizontal dashed line the "plateau value" 6r^^. Figure adapted from 11231 . 
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3.2 Effects of confinement on glassy polymers 



3.2.1 Brief overview of some experimental results 

During the past fifteen years the impact of geometric confinement on the glass tran- 
sition has received considerable attention. The progress in the field is described in 
several comprehensive reviews [23, 124i 11251 11261 l24l 11271 [1281 . The picture emerg- 
ing from these studies is that glass formers confined to nanoscopic dimensions may 
exhibit deviations from bulk behaviour due to the interplay of spatial restrictions and 
interfacial effects. The latter (nonuniversal) effects can result from particle-substrate 
interactions, confinement-induced changes of the liquid structure or polymer confor- 
mations, density variations, etc. llT3ll24l . and often appear to dominate the behaviour 
of the confined glass former II127L 

These interfacial effects have an important impact on the glass transition of thin 
polymer films. For films supported on a substrate, many studies 11291 11301 11311 11321 
[I33][l3l[l33[l36][l37][l38, m 141, 142 1431. though not all lfT44lfT45l. find 
reductions in Tg with decreasing film thickness (h) if the polymer-substrate attraction is 
weak. This is also the case for one of the most extensively studied systems, polystyrene 
on a variety of substrates [129, I3l\ \13^ [IMl ES] [136 137 139j HSl MB fT42l 
11431 . Here mesurements of surface relaxation after nanodeformation 014611 or of the 
positional dependence of Tg [139] further indicate that the depression of Tg is related to 
a "free-surface effect": Monomers near the free surface are expected to be more mobile 
because they feel less steric contraints than their peers in the bulk. This enhanced 
mobility should lead to reductions in Tg. 

Measurements of (the average) provide information on the dynamic reponse of 
the polymer films on the time scale associated with Tg which depends on experimental 
conditions, such as the cooling rate fl47|. Additionally, the full structural relaxation 
has also been explored by several techniques 1 126, 24J, in particular by dielectric spec- 
troscopy |[l22[l3i[T3l[l35][l48l[l4l][n2l[^ Dielectric spec- 
troscopy allows for a simultaneous measurement of Tg and the relaxation spectrum, 
even for thin films. A key finding of these studies is that the a process is broadened 
relative to the bulk. It is possible to interpret this broadening as a consequence of 
spatially heterogeneous dynamics in nanoconfinement. At the interfaces the segmental 
dynamics can be enhanced (e.g. at the free surface) or slowed down (e.g. at an attractive 
substrate) relative to the bulk. Chain segments in layers adjoining these interfacial lay- 
ers should also have their dynamics perturbed, albeit to a lesser extent, which will lead 
to a still weaker perturbation for the next layer, and so on. This interpretation implies 
that there is a smooth transition from interface-induced perturbations of the dynamics 
to bulk behaviour with increasing distance from the interfaces. Recently, this idea has 
been exploited to analyze dielectric spectra of nanostructured diblock copolymer melts 
M150I . and it is also consistent with the positional dependence of Tg found in Ml 391 . 

3.2.2 Simulation work on geometrically confined glass-forming polymer systems 

Simulation studies of confined glass formers support this view of an interface-induced 
gradient in relaxation ll32ll . Many of these studies utilize simple models, e.g. binary 
Uquids US] [Bll [155] [1561 or bead-spring polymer models 037] [BS] \159\ ]lM \lM\ 
[ISl[!S][lMl[IS][!Sl[!67][IM][!Ml[n3[ni][!^ (see however |[T75l 

11761 11771 11781 1179II for work on chemically realistic models), to explore the impact 
of confinement in thin films fT53\ \T55\ [BT] [Bi [153 [EQI [MD [Ml [Ml [Ml [Ml 
[M[M[M[l63[l70l[IlI][l72[l73[ni[n21[IlS[I^ pores HH or systems 
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Figure 9: Tc{h)/Tc and Tg{h)/To {T^ and Tg denote the bulk values) versus rescaled 
film thickness h/ho. MD results for supported films (open circles), free-standing films 
(shaded circles), and films confined between two smooth repulsive walls (filled circles) 
are compared to the glass transition temperatures Tg{h) of three studies: (i) Monte 
Carlo simulations of a lattice model for free-standing atactic polypropylene (PP) films 
01771 (crosses), (ii) Experiments of supported atactic PS films of low molecular weight 
(open squares) II132I . (iii) Experiments of supported, high-molecular weight PS films 
|fT29J (stars). The solid and dashed line show Eqs. Eq.|7]and Eq. |8] respectively. The 
vertical dotted line roughly indicates a fibn thickness of 10 nm. Figure adapted from 

mill. 



containing nanofillers II1801 11811 11561 . The simulations reveal a complex relaxation 
behaviour on approach to the glass transition of the confined glass former and also 
report shifts of Tg, qualitatively similar to the trends sketched above for experiments. 
In the following we illustrate these results by some examples. 

Figure |9]compares the reduction of Tg with film thickness found in experiments on 
supported PS films of low 111 321 and high molecular weight 111291 with the Tg shifts 
obtained from simulations of free-standing films 117011171111771 . supported films 1,1701 
11711 . and films confined between two substrates 01661 11671 . The simulations study a 
chemically realistic model of polypropylene 1 177] or flexible bead-spring models II 1661 
I167II1701[T71J . and use as a substrate completely smooth walls which are either purely 
repulsive or weakly attractive. Despite the obvious differences between experimental 
and computational systems — different polymers, absence or presence of a substrate, 
etc. — a master curve for the reduction of Tg can be constructed, if Tg{h) is scaled by 
the bulk value and h by a characteristic thickness ho that depends on the nature of 
the system, but only (very) weakly on molecular weight. This master curve allows 
us to compare the different systems, which is instructive in several respects. First, 
the close agreement of the PS data for low and high molecular weight suggests that 
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possible changes of the entanglement density II1371 1 1821 or chain conformation II183I 
11821 in thin films are probably not responsible for the Tg depression. Furthermore, 
the experimentally observed Tg shifts have stimulated several attempts to model this 
phenomenon theoretically lT32[lM[lMl[lM[l821[l88l[l89l[^ For 
instance, based on a percolation model for slowly relaxing domains Long and Lequeux 
II191I derive the formula (solid line in Fig.|9]l 



originally suggested by Keddie et al as an empmcal parametrization in their seminal 
study on supported PS films 11291 . An alternative parameterization was proposed by 
Kim et al 1(1561 (dashed Hne in Fig.|9]) 



and an attempt was undertaken to justify this formula by a viscoelastic capillary waves 
model 1 132, 184|. Figure |9] shows both formulas and reveals that, for a typical range 
of experimental film thicknesses, 15 < /j/Ziq ^ 300, it is hard to decide whether Eq.|7] 
or Eq. [8] is more accurate. However, including the simulation data for h <\Q nm, 
Eq. [8] appears to provide the better description of the Tg shift. Therefore, Eq.|8]will be 
employed later (Fig.fTSTi in an analysis of the local Tg of simulated polymer films. 

The explanation for the Tg shift in the simulations rests upon the impact that the 
boundaries exert, because the average behaviour of the film — and so its Tg — aggregates 
contributions from all layers in the film, and these layers have distinct properties. This 
can be illustrated by an analysis resolving structure and dynamics as a function of 
distance from the boundaries. For flexible bead-spring models Fig. [TO] shows two ex- 
amples of such an analysis: the layer-resolved incoherent scattering function (j)g{t,z) at 
q* (maximum of 5'(^); cf. Fig.fTTTi for a polymer melt surrounding a highly faceted, but 
nearly spherical filler particle with a structured and attractive surface (panel (a)) II180L 
and (p^t{t,z) for a polymer film supported by a smooth, weakly attractive substrate 
(panel (b)) flTOll. 

For the polymer melt surrounding the nanofiller the scattering function, averaged 
over all monomers (crosses in Fig. [TOt a)). displays features unfamiliar from the bulk. 
The a-relaxation seems to occur in two steps, as if there were two distinct processes, 
a fast one corresponding to a bulk-like phase far away from the filler and a slow one 
associated with interfacial relaxation. The layer-resolved analysis reveals that this in- 
terpretation is misleading. The strongly stretched tail of the average correlator results 
from the smooth gradient in the decay of 0|(f ,z) which slows down on approach to the 
filler particle. References [180 181] show that the amplitude of this tail can be tuned 
by the monomer-filler interaction. Strong attraction leads to a more pronounced tail; 
vanishing attraction suppresses the tail. In the latter case, the shape of the scattering 
function is bulk-like. 

Additional insight into the slow relaxation of particles in contact with a structured 
wall was obtained from studies of the self-part of the van Hove function M1551 11941 
1195 1. For times in the a regime and T ^ Tc, Gs{r,t) has a clear two-peak structure. The 
first peak reflects particles that remain trapped in their cages, while the second peak, 
located on a length scale corresponding to the wall structure, reveals that particles have 
migrated to neighbouring wells. Evidence for this kind of "hopping motion" is found 
for a binary LJ mixture 1, 1551 . for a polymer melt adsorbed on a structured surface 
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Figure 10: (a) Layer-resolved incoherent scattering function (j)y{t,z) at T — 0.4 and 
^ = 7.08 (= maximum of 5'(^)) for a bead-spring polymer melt surrounding an icosa- 
hedral filler particle. The chain length is N = 20, the melt density is p = 1 (dashed 
horizontal line in the inset), and the filler attracts the monomers more strongly than 
they attract each other in the bulk. The solid lines and the squares show (pg{t,z) for 
different distances z from the surface of the filler particle. The location of the first 
two layers is illustrated in the inset which depicts the monomer density profile p {z/Rg) 
(Rg ~ 2.17). The crosses indicate the average over all layers. Figure adapted from ll32ll . 
(b) ^g{t,z) at r = 0.42 and q = 6.9 (w maximum of S{q)) for a bead-spring polymer 
melt in a supported film of thickness h = 20.3 (Tc ~ 0.392). z denotes the distance from 
the (left) wall. 0^(f ,z) is obtained as an a^grage over all monomers which remain for 
all times shown in a layer of width Az ~ 2 and centered at z- The average behaviour 
of the film (average over all layers) is indicated by crosses and the bulk data by filled 
circles. Inset: Monomer density profile p{z) versus z- The layers for which <^|(f,z) is 
shown in the main figure, are indicated. Figure adapted from II170I . 
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Figure 11: Static structure factor S{q) for a bead-spring melt (A^ — 10) in the bulk 
(circles) and in supported film geometry (lines) at T = 0.44 {Tq w 0.405 in the bulk; 
Tc{h ~ 7) « 0.361 in the film [1701). For the film S{q) is shown respectively for layers 
in film center, at wall and at the free surface. At the free surface the steep rise of 
S{q) for small moduli q of the wave vector is due to capillary waves. For the bulk the 
compressibility plateau k^TpKj is indicated by a horizontal dotted line (Kj denotes 
the isothermal compressibility). Inset; Same data as in the main figure, but in a log-log 
representation. The small-g' behaviour expected from capillary wave theory is depicted 
by a dashed line (7 is the surface tension). 

111941 or for model of liquid toluene confined in cylindrical mesopores [195]. This 
relaxation mechanism could be generic when a liquid may lock into registry with the 
surface topography, leading to a mechanism of structural slowing down which coexists 
with glassy arrest at low T II 1961 [197 J 

A similar locking is not possible at the free surface or a smooth, repulsive or weakly 
attractive wall. In thin films, one therefore expects enhanced dynamics relative to the 
bulk. Figure [TOl'b) shows that this expectation is borne out for flexible bead-spring 
models. For these models the free and smooth interfaces also create enviroments for 
nearby monomers which tend to reduce the cage effect: Local spatial correlations on 
the scale q* of the maximum of 5'(^) are weaker in the films than in the bulk at the 
same temperature (Fig. fTTT i. This is an important contributing factor to the enhanced 
monomer dynamics found at both interfaces. 

Figure [TOl'b) also shows that the interface-induced deviations from bulk dynamics 
continuously turn into bulk-like relaxation with increasing distance from the bound- 
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Figure 12: Temperature dependence of the layer-resolved a relaxation time Ta(z, T) in 
a supported polymer film (A^ = 10; 7^ « 0.405 in the bulk). The simulated system is 
the same as in Fig. fTOl b). The solid lines at the free surface show Eq.[TT](no adjustable 
parameter). The horizontal dashed lines indicate the bulk values for Ta at the respective 
T . Figure adapted from 01711 . 
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aries. The range of this crossover grows on cooHng. Figure [12] illustrates this point 
by an analysis of the temperature dependence of the layer-resolved a relaxation time, 
"^aiZjT). At high temperature Ta{z,T) slightly deviates from the bulk value near the 
interfaces. With decreasing T the interface-induced enhancement of the dynamics in- 
creasingly penetrates into the film and eventually propagates across the entire system 
for sufficiently low T. This spatial dependence of Ta{z,T) can be modeled if two as- 
sumptions are made: (i) The average Tc{h) of the film is given by Eq. [8]and can be 
written as a "democratic average" 11931 of the local Tdz). That is. 



which gives 

Ul+ho/z) 

"^'^'^ = {l+ho/2zV • ^^'^ 

(ii) The second assumption is that the sole effect of the interface is to shift Tc from the 
bulk value to Tc{z), whereas all other parameters determining the T dependence of Ta 
remain the same as in the bulk. Here we model this T dependence by the MCT power 
law 

Figure [T2l demonstrates that Eq. [TT] yields a reasonable description of the increase of 
Ta from the free surface to the center of the film (solid lines in the figure) |T71 1. The 
decrease of the a relaxation time on approach to the free surface therefore corresponds 
to a decrease of the local glass transition temperature, in qualitative agreement with the 
experimental results of 1 139 1 

The propagation of enhanced or reduced mobility from the boundary toward the 
interior of the film has also been observed in other simulations on freely-standing 
[T58| and supported polymer films 11641 . These studies carried out a cluster analy- 
sis, of highly mobile monomers in the case of the freely-standing film and of immobile 
monomers for the supported film. In both cases, it was found that clusters start at the 
interface and penetrate into the film. 

In summary, simulation studies suggest that confined (polymeric) liquids display 
complex relaxation behaviour on approach to the glass transition because of the in- 
terplay of bulk-like slowing down of the dynamics and interfacial effects. Interfaces 
can enhance or retard the relaxation relative to the bulk. Enhanced relaxation may be 
expected for smooth or free interfaces, whereas strong particle-substrate attraction or 
particle caging in cavities of the substrate tend to slow down the dynamics. These 
interface-induced perturbations smoothly transition from the boundaries to the interior 
of the confined liquid. The range of this gradient grows on cooling so that the perbur- 
bations can propagate across the entire liquid for sufficiently strong confinement or low 
T. Similar interfacial effects may also be important for the analysis of other problems, 
for instance, for solvent evaporation from (spincoated) polymer films II173II1741 or for 
the modeling of the hydrodynamic boundary conditions in microfluidic devices [1981. 
Therefore, it appears that those theoretical approaches, which treat the interplay of 
boundary effects and glassy slowing-down of the dynamics on the same (microscopic) 
footing, are most likely to advance our understanding in this field. Folding in boundary 
effects is certainly a major challenge [199. 196. 1971 . Simulations of model systems — 
as those reported here — should be helpful for the development of such theories. 
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4 Studies of the glassy state 



Molecular simulations are intrinsically limited in terms of time scales. Therefore study- 
ing the glassy state, in which the relaxation time scales are by definition extremely 
large, could seem to be out of reach for this kind of numerical approach. Paradoxically, 
however, the fact that the intrinsic time scales of the system are large brings the simu- 
lations very close to actual experiments. The important fact is that, deep in the glassy 
state, segmental relaxation times |200| greatly exceed both the experimental and the 
simulation time scales. For both laboratory and computer experiments, a time window 
is probed where segmental polymer motions are frozen and hence non-equilibrium phe- 
nomena associated with the glassy state are observed. The main issue is then rather a 
matter of sample preparation, which is usually done with much faster quenching rates 
in computer experiments (cf. Sect. 13.1. 11 1. As we will see below, it turns out that the 
mechanical behaviour is not strongly affected, at least at a qualitative level, and that 
simulations can therefore be used to analyse the structure and mechanical response of 
the glassy state with some confidence. 

4.1 Small deformations 

We start our discussion by considering the small strain, elastic part of the mechanical 
response of a polymer glass, which can be described at the macroscopic scale by linear 
elasticity. It is now well recognized that, in spite of their homogeneity in density, glassy 
materials are heterogeneous at the nanoscale in terms of their elastic properties. This 
heterogeneity is reflected indirectly in some vibrational properties (e.g., the Boson peak 
112011 '). and is most obviously revealed by simulation, that allows one to compute elas- 
tic constants and study elastic response at various scales [202 j . For polymer glasses, 
the first study of that kind was presented in refs M2031 12041 . following earlier studies 
that indicated size dependent results for the elastic constants of nanometric systems 
ll2l)5l . This study, which uses a simple definition of local elastic constants based on 
a local calculation of the usual fluctuation formulae 02O41 12061 . shows clearly that the 
material is inhomogeneous at scales of the order of 5 to 10 monomeric sizes. Small 
regions displaying negative elastic constants and therefore would be unstable if they 
were not surrounded by regions with large local moduli. These aspects are not specific 
to polymer systems, but can be observed in simple molecular or metallic glass formers 
II207II2O2J . and using more sophisticated definitions of the local elastic moduli, which 
points to their universal character. However, they are particularly important when the 
polymer glass is cast in the form of a nanostructure, as the mechanical failure of such 
structures will be strongly affected by the presence of elastic heterogeneities. 

While the response of a polymer glass for deformations smaller than a few percent 
can be described as elastic from a macroscopic viewpoint, a detailed microscopic stud- 
ies show that irreversible rearrangements take place at low temperature even for very 
small deformations 12081 . This irreversible behaviour is responsible for a dissipative 
part or the mechanical response, sometimes described as anelasticity. Ultimately, this 
dissipation is associated with strongly localized plastic events that involve only a few 
monomers, and can be thought as the precursors of the dissipative processes that take 
place under larger deformation. 

Two specificities of polymer glasses, already mentioned in Sect. 13.2.21 are the in- 
terest for thin films and the ability to introduce various types of additives, from small 
molecules to nanoparticle fillers. Relatively few studies are concerned with the thin 
film modifications of the glassy state itself. Jain and de Pablo I209J showed that the 
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vibrational density of states is modified compared to the same glass in the bulk. In 
particular, the usual excess of low frequency modes observed in glasses (the so called 
boson peak) is enhanced in free standing thin films. This increase of the soft modes, low 
frequency part of the spectrum is consistent with a negative shift of the glass transition 
temperature. As the strength of the boson peak is in general, connected to a decreased 
fragility (in the sense of Angell's classification of glasses, see |12|, this points to a 
reduction in fragility with film thickness. 

Many studies on the other hand have been devoted to the modification of the glassy 
state under the influence of nanometric filler particles Ill80lll8l1l210l . As discussed in 
Sect. 13.2.21 (see Fig. [TOl i. the vicinity of fillers was shown to modify the dynamics of 
the polymer, with a slowing down and a shift of the glass transition temperature in the 
case of attractive interactions. From the mechanical point of view, it was also shown 
that the presence of attractive fillers results in the presence of a layer with enhanced 
mechanical properties (higher local moduli) 1121 1! 121 01 and modifies the structure of 
the entanglement network II212I . The spatial extension of this layer is of the order 
of a few monomer diameters and does not depend on the particle radius as soon as 
the latter exceeds 2 to 3 monomer size. As a result the mechanical properties of the 
composite are enhanced, however for a well dispersed composite this enhancement 
remains moderate. Note that the existence of a glassy layer around the nanoparticles is 
expected to affect the properties of the nanocomposite mostly in the region of the glass 
transition, where the contrast between this layer and the rest of the polymer is maximal 
II2131 12141 . As for polymer films, the fragility (inferred from boson peak strength and 
position) appears to be decreased with respect to the pure polymer 

Further studies of weak deformation include the influence of an external strain on 
the microscopic dynamics |215| (in fact the latter study also extends to large, uniax- 
ial deformations) , and the evolution of the creep compliance upon aging ll2T6j . It is 
found a finite strain rate accelerates the segmental dynamics -measured by the bond 
orientation relaxation time- , whether the corresponding strain is positive (dilation) or 
negative (compression). This symmetry between compression and dilation is indica- 
tive of a stress induced modification of the potential energy landscape, with the applied 
stress lowering the local barriers and accelerating the caging dynamics. Free volume 
considerations, on the other hand, would not explain such a symmetry f2T5\. Acceler- 
ation of the dynamics under load (sometimes described as "mechanical rejuvenation") 
is also observed in the aging study of ref; [216|, in which the creep compliance J{t,tn.) 
of a simulated polymer is also shown to be well accounted for by the classical descrip- 
tion of Struik 02171 . in which the creep compliance is described as a scaling function 

J{tgff/tn,) where f,j. is the aging time, and f^yy — Jq i^'T+p'^ This result again 

shows that, in spite of the wide difference in time scales, the phenomenology of glassy 
polymers is reproduced by simulation work. It should also be noticed 1,218] that the 
relation between segmental mobility and strain rate is also dependent on the general 
strain history, so that, in contrast to the usual assumptions of Eyring's theory, no sim- 
ple mechanical variable can relate to the local mobility. Recent theoretical approaches 
mOOI . based on nonlinear Langevin equation description of segmental dynamics, go 
beyond the simple Eyring description and account for this aging behaviour 

Finally, a rich and promising domain for simulation is the study of the complex 
effects of small molecules on the glass transition and on the glassy state. Such small 
molecules can be described as solvent, plasticizers or antiplasticizers, and in general 
their presence leads to an acceleration of the dynamics and a decrease in the glass 
transition More precisely, plasticization of the polymer means that the solvent does not 
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only decrease Tg, but also softens the polymer glass by reducing the elastic moduli. 
Besides this normally encountered case, there are also systems in which the solvent 
antiplasticizes the polymer. That is, the solvent decreases Tg, but increases the elastic 
moduli of the polymer glass. The origin of this antiplasticizing effect has recently been 
studied by simulations |219, 220 1, and appears to be due to a more efficient packing 
associated with solvent molecules smaller than the monomers. Small molecule fit into 
the "holes" of the polymer melt, and increase the elastic stiffness. However, their 
mobility also facilitates segmental mobility. Both for antiplasticizers and plasticizers, 
a strongly heterogeneous dynamics of the solvent molecules has been reported 12201 
11231 . The corresponding dynamical correlation length is however decreased in the 
presence of antiplasticizres, with a fragility of the system that is also decreased. 

4.2 Large deformations and strain hardening 

The mechanical properties of glassy polymers at large deformations were first inves- 
tigated in the pioneering work of Rottler and Robbins [221 , 222 223 224], and their 
results were subsequently reproduced in a number of studies Ii76. i225|. Beyond the 
peak stress, whose actual amplitude depends on the strain rate in a logarithmic way, 
as in simple glasses, the stress strain behaviour depends on the type of sollicitation. 
In pure shear or in simulations of a tensile test under triaxial conditions, the plastic 
flow proceeds through a cavitation of cavities and subsequent formation of fibrils. The 
drawing of these fibrils is essentially at constant stress. This phenomenon is discussed 
in detail by Rottler |27 1 Under uniaxial, almost volume conserving conditions (Poisson 
ratio close to 1/2), a marked strain hardening is observed, i.e. the plastic stress increases 
with strain with a dependency that is close to linear. The corresponding deformation, 
however, is very far from elastic, and essentially non recoverable. 

This peculiar property of polymer glasses, which strongly contributes to their prac- 
tical applications, is the existence of the strain hardening regime at large deformation. 
Large strain-hardening effect prevents the strain localization and leads to a tough re- 
sponse of the polymeric material; the material breaks only after a significant plastic 
strain, as in the case for polycarbonate. For brittle polymers, as atactic polystyrene, 
the strain-hardening effect is usually too weak, which leads to a brittle fracture due to 
crazes formation within a few percent of strain. Understanding the microscopic ori- 
gin of the strain-hardening effect provides a strategy for new tailor-made polymeric 
materials. 

The microscopic origin of this regime has been clarified only recently, thanks - to a 
large extent - to numerical simulations. Both the generic aspects of this phenomenon, 
and the way strain-hard and strain-soft polymer materials differ in terms of segmental 
mobility and energetics, have been investigated by extensive MD simulations. 

Before we discuss simulation results, let us recall the popular rubber elasticity mod- 
els [226 1 which predicts that the hardening modulus G/, varies linearly with tempera- 
ture T and entanglement density Pe, G/, = pek^T. Experimental values, however, are 
two orders of magnitude larger, and, more important, show decrease with increasing 
T 02271 12281 . To take the dissipative nature of the plastic deformation into account 
the rubber elasticity model describes the stress-strain relation in the strain-hardening 
regime with the help of two contributions 02261 . 

CT = C7j' + C2(A2-A-1) (12) 

The first part is a constant dissipative stress Gy , due to the presence of energy 
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Figure 13: Snapshots of three configurations of a polymer glass: initial undeformed 
state (up), uniaxially deformed configuration (left), and triaxially deformed configura- 
tion (right). Both deformed configurations correspond to a strain of 50% along the z 
axis of the simulation box. The simulation results are obtained from a flexible bead- 
spring model (A^ = 100) at T = 0.2 (T^g w 0.43). Figure taken from |i76|. 
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barriers. The second is the strain-dependent part, which is thought to be described by 
rubber-elasticity theory and represents the strain-hardening effect. In this description 
the strain-hardening modulus is not affected by thermally activated processes. 

However, recent experiments demonstrate that such a description for the strain- 
hardening part is invalid; the strain-hardening modulus has characteristics of a thermally- 
activated process, and decreases for higher temperatures II229II . Secondly, at higher 
strain rate the strain-hardening modulus increases II230I . although for some polymers 
the dependency on strain rate is rather weak [231]. This increase can be interpreted in a 
very elementary manner within a barrier-crossing, Eyring like picture, which however 
as discussed above is known to be oversimplified. Finally, the external pressure affects 
the strain-hardening modulus as well. A higher external pressure leads to an increase 
in the strain-hardening modulus 12321 . Again, this behaviour is typical for thermally- 
activated processes. All these three observations on the strain-hardening modulus are 
not present within the classical rubber theory. The failure of rubber-elasticity theory 
is due to the essential difference between a rubbery state, with many possible chain 
conformations between cross-links, and a non-ergodic glassy state where chain confor- 
mations are practically frozen, and transitions between different conformations are not 
possible. Deformation of the polymer glass facilitates these transitions and is accom- 
panied with the energy dissipation. The work for this irreversible energy dissipation is 
reflected in a dissipative stress, which is absent in the rubber-elasticity theory. The dis- 
sipative nature of a polymer strain-hardening is confirmed by experiments. It is found 
that for PC and PS at large (> 15 - 30%) strains Il233ll20ll234ll ) more work is dissipated 
through heat than converted into internal energy. 

The increase in stress for larger strains in combination with the dissipative nature 
of the stress in the strain-hardening regime suggest that there is an increase in the rate 
of energy dissipation, i.e., more energy per unit of strain is needed for more stretched 
samples to stretch them further This picture is supported by the computer simula- 
tions of Hoy and Robbins |235|, who showed that the dissipative stress increases with 
larger strain and that at zero temperature the stress was directly correlated to the rate of 
changes in Lennard- Jones (LJ) binding. Their more recent simulations of polymer toy 
models also demonstrated that most of the stress at large strains is due to dissipation. 
I1235J 12361. 

If the rubber elasticity theory is invalid, what is then the polymer-specific part of 
the strain-hardening modulus? Simulations of the mechanical deformation of simple 
molecular glasses show no strain hardening |237|. Hence, for short polymer chains 
the amount of strain hardening is expected to be small as well. Moreover, experiments 
12281 and simulations |238| show that the strain hardening modulus is positively cor- 
related with the entanglement density. Therefore one expects the strain hardening phe- 
nomenon to disappear progressively as the chain length falls below the entanglement 
threshold. Simulations show that this is indeed the case. In a model polymer glass. Hoy 
and Robbins [235] demonstrated that there is a gradual increase in the strain-hardening 
modulus as a function of chain length for chains up to about the entanglement length. 
For longer chains saturation in the modulus occurs. Similarly, in MD simulations of 
atactic polystyrene 1,48] only a weak strain hardening was observed for chain lengths of 
80 monomers, below the experimentally observed entanglement length of about 128- 
139 monomers ||239ll240l . 

In the simulations of ref. | 235| it was also observed that the strain hardening is 
more correlated with the change in the end-to-end distance of the polymer chains than 
with the change in the global sample size. If the sample size is decreasing while the 
end-to-end distance does not, then the stress does not increase. 
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After the fundamental aspects of strain hardening have been clarified using simula- 
tions of simple, coarse-grained models, more detailed studies using realistic models are 
required to understand the influence of chain architecture. For example, it is known that 
polymers with a larger persistence length often have a higher strain-hardening modulus 
11226 1 . It was even shown by molecular-dynamics (MD) simulations that if the persis- 
tence length of a polyethylene-like polymer is artificially increased by changing the 
trans -to-gauche ratio, the strain-hardening modulus of the resulting material increases 
as well [2411. That example illustrates that the strain-hardening modulus depends on 
the conformation of the chain, which is frozen in the glassy state. The changes in this 
conformation are associated with plastic events, that involve a collective nanoscale seg- 
mental dynamics of individual or neighbouring chains. Understanding the differences 
in local mobility and dynamics for chemically different polymers under large defor- 
mations, and how it relates to the different mechanical characteristics, is therefore an 
important question, which has been studied in relatively few cases. 

As an example of such studies, and of the possibilities offered by simulation, we 
consider the extensive studies in refs. M242i l48l 12321 1243 1 on two glassy polymers that 
vary greatly in their strain-hardening moduli, viz. polystyrene (PS) and polycarbonate 
(PC, of which the modulus is more than a factor of two higher |244|). Molecular- 
dynamics simulations Il242ll48l[2321 have reproduced these experimental findings qual- 
itatively, with a strain-hardening modulus of polystyrene that is much lower than that of 
polycarbonate, as shown in figure [141 They also allow one to identify the microscopic 
mechanisms responsible for the observed difference in the strain-hardening modulus. 
In particular, they show the rate of non-affine displacements (or local plastic events) 
increases with larger strain, and that the increase is larger for PC. The non-affine dis- 
placements of particles are due to restrictions and hindrances, in particular from cova- 
lent and steric interactions. These restrictions are extremely important at the scale of 
the covalent bond. If particles would displace affinely, then the equilibrium value of 
this chemical bond would be excessively disturbed. To circumvent the bond stretch, 
the bond vector will not move affinely with the deformation. At much larger length 
scales the situation is different. For a long chain the internal conformation can be ad- 
justed, while still obeying to the equilibrium length of the covalent bond. The effective 
spring constant associated with the end-to-end distance is much less stiff. Moreover, 
due to the glassy state the relaxation time of spontaneous rearrangements at the scale of 
the whole chain greatly exceeds experimental and simulation time scales . Therefore, 
the end-to-end distance cannot adjust back towards the equilibrium value by means of 
spontaneous relaxations. Hence the end-to-end distance will follow the imposed defor- 
mation much more affinely. So we expect a more affine response, as we probe larger 
length scales. 

This expectation is borne out by a detailed study of the affine character of the de- 
formation of the polymer chain, as a function of the internal distance along the chain 
||243J . At the scale of 100 chemical bonds, it is found that the deformation is essen- 
tially affine up to about 15% strain for both PS and PC, while it is strongly non affine 
already at this strain at the scale of 30 bonds, see Fig. [15] For larger strains the relative 
deviation from affinity becomes larger for even long enough (100 backbone bonds) 
chains. This effect is present in both PS and PC, the difference being in the magni- 
tude of the effect. For an internal distance of 30 monomeric units, the deformation 
of the polycarbonate chain is only 30% of the affine value at a strain of 50%, while 
in polystyrene under the same conditions the chain deformation is 50% of the affine 
value. This is because at the scale of the Kuhn length a chain cannot be stretched any 
further As the Kuhn length of polycarbonate is already larger than that of PS and the 
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Figure 14: The von Mises equivalent true stress a vs. strain Ceng for PS and PC. Solid 
lines are fits to Eq.[T2l Fit range is Ceng = 0.3-0.8. Note that the strain-hardening mod- 
ulus Gh for PC (19 MPa) is almost twice that of PS (1 1 MPa), while their extrapolated 
offset yield values Oy are about the same (PC: 88 MPa; PS: 86 MPa). Adapted from 
ref. ESI 
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total non-affine displacement of PC is to a large extent determined by the backbone 
atoms, the increase in the effective stiffness length leads to an increase in non-affine 
displacement, more energy dissipation and hence a higher strain-hardening effect. For 
polystyrene the Kuhn length is small and the major part of non-affine displacement 
is not caused by the backbone. Hence the expected increase in the effective stiffness 
length during deformation does not lead to a substantial increase in plastic flow, so that 
at moderate strains PS behaves more like a simple glass without strain hardening, as 
opposed to polycarbonate. 

5 Perspectives 

The present review has attempted to describe some aspects of the recent progresses in 
our understanding of glassy polymers that have been obtained from simulation work, 
with a particular focus on microscopic dynamic, glass transition and mechanical prop- 
erties. The results obtained in the past ten years show that the simulations have reached 
a state of maturity that allows them to address and clarify important issues in the field, 
in spite of the obvious limitations in terms of length and time scales. While the il- 
lustrative examples we have chosen were relatively simple, they show that the general 
phenomenology of the glass transition and the properties of the glassy state can be 
accounted for using molecular simulations, in spite of the relatively small time scale 
that can be considered in such simulations. We would like to close the discussion by 
suggesting some directions for future research, which we believe will be soon - or 
are already - within the capabilities of state of the art simulations, and correspond to 
problems of practical or fundamental interest. 

A major ingredient of the physics of glassy polymers - and more generally of glassy 
systems - is the so called "time temperature superposition principle" which assumes 
that relaxation processes can be rescaled by using a single, temperature dependant re- 
laxation time Ta(r). While this "principle" is widely used to produce frequency depen- 
dent relaxation function from data taken at various temperature, it is also well known to 
be applicable at best in the frequency region corresponding to the a peak, presumably 
only above Tq. It is a challenge to understand the molecular mechanisms underlying 
this superposition principle and the deviations thereof 01 181 . A prominent feature of 
polymeric glasses is the observation of secondary relaxations at higher frequency, the 
Johari-Goldstein j3 processes (see Sect. l3.1.2] i. The associated relaxation times merge 
with the a relaxation times typically at temperatures of O.STg, and the associated time 
scales are of the order of 10^^-10^^^ at such temperatures. Such time scales are be- 
ginning to be in the range attainable by MD simulations, a fact that should allow one 
to clarify the molecular origins of secondary relaxations, their relation to the existence 
of side groups, to torsional modes, and other interpretations that have been put forward 
02451 . Molecular simulation, by the flexibility it offers to block specific motions and 
to control chain properties at various scales of coarse graining, should also allow one 
to investigate the putative relations between these secondary relaxations, fragility, and 
mechanical properties. 

Nanocomposites and thin films have already been the subject of a number of simu- 
lation studies (see Sect. l3.2.2b . Still, the precise mechanisms that give rise to reinforce- 
ment and nonlinear behaviour in nanocomposites, and to speed up of the dynamics in 
thin films, remain to be elucidated in detail. In particular, for entangled polymers it 
can be expected that the presence of interfaces modifies both the characteristics of the 
entanglement network and the monomeric friction coefficients. The interplay between 
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these different aspects could, in principle, be elucidated in simulations of entangled 
systems. Again, the simulation times needed to explore such phenomena are extremely 
large, but, with the development of advanced simulation algorithms — as those alluded 
to in Sect. |2] — and the increase of computer power they are starting to be within the 
attainable range. 

Finally, we note that many polymer materials of practical interest display a mixed 
structure, being either partially crystalline with an amorphous fraction, or chemically 
inhomogeneous with phases (or microphases in the case of block copolymers) display- 
ing different mechanical properties, e.g. glassy and rubbery. The properties of the 
resulting nanocomposites have, up to now, received Uttle attention from the standpoint 
of molecular simulation (see however F1 08II1071I246II ). They are likely to be strongly 
influenced, especially for large deformations, by the interplay between chain architec- 
ture and material nanostructure. Another interesting phenomenon of that kind is strain 
induced crystallisation, which plays an important role in the strain hardening of some 
specific polymer. The capability of molecular simulation to describe polymer crystal- 
lization has been recently demonstrated ITllSl 12471 . so that the study of strain induced 
crystallisation -as well as studies of semi crystalline phases appears to be an interesting 
goal for the future. 
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